Elliptic flow of colored glass in high energy heavy ion collisions 
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We compute the elliptic flow generated by classical gluon fields in a high energy nuclear collision. 
The classical gluon fields are described by a typical momentum scale, the saturation scale A s , which 
is, for RHIC energies, of the order of 1 — 2 GeV. A significant elliptic flow is generated only over 
time scales on the order of the system size R. The flow is dominated by soft modes pt ~ A s /4 
which linearize at very late times r~fl> 1/A 3 . We discuss the implications of our result for the 
theoretical interpretation of the RHIC data. 



INTRODUCTION 

The collective flow of excited nuclear matter has been 
an important tool in attempts to extract the nuclear 
equation of state ever since the early days of heavy 
ion collision experiments Q. Measurements of collec- 
tive flow at the Relativistic Heavy Ion Collider (RHIC) 
may provide insight into the excited partonic matter, of- 
ten called a Quark Gluon Plasma (QGP), produced in 
high energy heavy ion collisions j^]. In particular, the 
azimuthal anisotropy in the transverse momentum dis- 
tribution has been proposed as a sensitive probe of the 
hot and dense matter produced in ultra-relativistic heavy 
ion collisions || • A measure of the azimuthal anisotropy 
is the second Fourier coefficient of the azimuthal distri- 
bution, the elliptic flow parameter V2- Its definition is B 
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The elliptic flow for non-central collisions is believed to 
be sensitive to the early evolution of the system S . 

The first measurements of elliptic flow from RHIC, at 
center of mass energy ^snn — 130 GeV, have been re- 
ported recently ||. Hydrodynamic model calculations 
provide good agreement, for large centralities, and for 
particular initial conditions and equations of state, with 
the measured centrality dependence of the data || . The 
agreement at smaller centralities is less good, perhaps re- 
flecting the breakdown of a hydrodynamic description in 
smaller systems ||. Hydrodynamic models also agree 
well with the pr dependence of the unintegrated (see 
Eq. (jj])) elliptic flow parameter V2(pt) up to 1.5 GeV/c 
at mid-rapidity [0. However, above 1.5 GeV, the exper- 
imental distribution appears to saturate, while the hy- 
drodynamic model distributions continue to rise j^]. Jet 
quenching scenarios to explain this saturated behavior of 
V2(pt) at large pt |l0| appear to disagree quantitatively 
with the data. Partonic transport models including only 
elastic gluon-gluon scattering require large cross sections 
or large initial gluon number to obtain significant elliptic 
flow [[l). 



In this letter, we compute the contribution to vi at 
central rapidities from the strong color fields generated 
in the initial instants after the heavy ion collision. These 
are generated as follows. At high energies, or equiva- 
lently, at small Broken x, the parton density in a nucleus 
grows very rapidly and saturates eventually ]l^ | form- 
ing a Color Glass Condensate jl| (CGC). The CGC is 
characterized by the color charge squared per unit area 
A 2 , which grows with energy, centrality and the size of the 
nuclei. Estimates for RHIC give A s ~ 1 - 2 GeV g|. 
For a recent review of the CGC model and additional 
references, see Ref. (Q. Since the occupation number of 
gluons in the CGC is large, / ~ l/as(A 2 ) > 1, classi- 
cal methods can be applied to study gluon production in 
heavy ion collisions at high energies [jl8| [r]] . The energy 
and number |2^, |lj of gluons produced were computed 
numerically for an SU(2) Yang-Mills gauge theory and 
recently extended to the SU(3) case |p2| . We have con- 
firmed that strong electric and magnetic fields of order 
1/ 'as are generated in a time r ~ 1/A S after the collision. 

The classical Yang-Mills approach may be applied to 
compute elliptic flow in a nuclear collision 0]. For pe- 
ripheral nuclear collisions, the interaction region is a two- 
dimensional almond shaped region, with the x axis lying 
along the impact parameter axis and the y direction per- 
pendicular to it and to the beam direction. We will show 
that even though large electric and magnetic fields (and 
the corresponding transverse components of the pressure 
in the x and y directions) are generated over very short 
time scales r ~ 1/A S , the significant differences in the 
pressures, responsible for elliptic flow, are only built up 
over much longer time scales r ~ R. Moreover, the ellip- 
tic flow is generated by soft modes pt ~ A s /4. Our result 
has important consequences for the theoretical interpre- 
tation of the RHIC data-these will be discussed later in 
the text. 



NUMERICAL METHOD 

We now discuss our numerical computation of ellip- 
tic flow. As in our earlier work, we assume strict boost 
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invariance. The dynamics is then that of a Yang-Mills 
gauge field coupled to an adjoint scalar field in 2+1- 
dimensions. For a numerical solution we use lattice dis- 
cretization. The discrctized theory is described by a 
Kogut-Susskind Hamiltonian pp| . 

In previous work, we studied gluon production in cen- 
tral collisions of very large nuclei and therefore assumed a 
uniform color charge distribution (A s = constant) in the 
transverse plane. To study effects of anisotropy and spa- 
tial inhomogeneity, we shall consider a finite nucleus. Wc 
shall impose suitable neutrality conditions on the distri- 
bution of color sources p3 to prevent gluon production 
at large distances outside the nucleus. 

To this end, we model a nucleus as a sphere of radius 
R 7 filled with randomly distributed nucleons of radius 
I w 1 fm. For a gold nucleus, R « 6.5 fm. The color 
charge distribution within a nucleon is generated as fol- 
lows. First, we generate (throughout the transverse plane 
of a nucleon) a random uncorrelated Gaussian distribu- 
tion p a (r) (a being the adjoint color index and r the 
transverse-plane position vector) , obeying the relation 

where the () average is over the ensemble of nucleons. 
Next, we remove the monopole and dipole components of 
the distribution by superimposing the distribution with 
the appropriate homogeneous contribution; first of the 
color charge, then of the color dipole moment. For a suffi- 
ciently fine lattice discretization, this procedure does not 
result in a significant change in the average magnitude of 
the random charge distribution. Since the color charges 
of different nucleons are uncorrelated, the resulting nu- 
clear color charge squared per unit area has a position- 
dependent magnitude, 
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, where r is the transverse radial coordinate relative to 
the beam axis through the center of the nucleus and I 
is the nucleon diameter. We can adjust A n to ensure 
the central nuclear color charge squared per unit area 
A 2 = Ag(0) has a desired value. 

Once the color charge distributions of the incoming 
nuclei are determined, the corresponding classical gauge 
fields can be computed. The initial conditions for the 
gauge fields in the overlap region between the nuclei are 
obtained as discussed previously p0| . For each config- 
uration of color charges sampled, Hamilton's equations 
are solved on the lattice for the gauge fields and their 
canonical momenta as a function of the proper time r. 



RESULTS 

We first compute the momentum anisotropy parame- 
ter a (defined in Fig. Q|) as a function of the proper time 
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FIG. 1: The momentum anisotropy parameter a — (T xx — 
T yy )/{T XX + T yy ) for a peripheral nuclear collision corre- 
sponding to impact parameter b/2R = 0.5 is plotted versus 
the proper time r in units of the nuclear radius R for two 
different values of A s oR. 



t p5j . The results for values of the external parame- 
ter A s0 i? = 18.5 and A s0 i?. = 74 (spanning the RHIC- 
LHC range of energies) are shown in Fig. [|. We observe 
that a rises gradually saturating at a ~ 1% at a proper 
time on the order of the size of the system. The time 
required to develop an anisotropy is clearly much larger 
than the characteristic time ~ 1 / A s o associated with non- 
linearities in the system. 

The calculation of elliptic flow, defined by Eq. ([!]), in- 
volves determining the gluon number, a quantity whose 
meaning is ambiguous outside a free theory. Closely fol- 

, we resolve this ambiguity 



lowing our earlier work |21 



by computing the number in two different ways; directly 
in Coulomb Gauge (CG) and by solving a system of relax- 
ation (cooling) equations for the fields. Both definitions 
give the usual particle number in the case of a free theory. 
We expect the two to be in good agreement for a weakly 
coupled theory. Wherever the two disagree strongly, we 
should not trust either. 

It is easy to show that w 2 N, N being the total gluon 
number, can be reconstructed from the cooling time his- 
tory of T xx — T yy , just as N can be reconstructed from 
that of the energy functional |2l] : 



v 2 N 
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(T xx (t)-T yv {t)). 



(2) 



This expression for v-i N is manifestly gauge invariant. 

In contradistinction to the gluon number, an estimate 
of V2 involves both the fields and their conjugate mo- 



menta. Indeed, consider the expression for T x 
our system: 



T 

1 it- 



T xx - T yy = J d 2 x ± [E 2 y -E 2 X + (D x $f - [Dy^f] , 
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where E is the chromo-electric field, $ the adjoint scalar 
field, and Di the covariant derivative, In the weak- 
coupling limit Di reduces to 9^, the usual derivative. In 
that limit, the first two terms in T xx — T yy only involve 
the conjugate momenta of gluons polarized in the trans- 
verse plane, while the last two terms only depend on 
the fields of gluons whose polarization is perpendicular 
to that plane. Since it is not a priori obvious that the 
two polarizations contribute equally to i>2 , both the fields 
and the conjugate momenta should be computed. For the 
cooling method, relaxation equations for conjugate mo- 
menta require that the usual relation between the mo- 
menta and proper time derivatives of the fields hold at 
all cooling times . 

In Figure |^, we compare, for a fixed impact parameter 
(b/2R = 0.5), the values of «2 obtained by the different 
methods. In the cooling approach, i>2 can be computed 
first by considering only the potential part of T xx — T yy 
in Eq. ^| and then assuming an equal contribution from 
the kinetic part. As seen in Figure such an equality 
does not hold until very late times. There is a significant 
difference at early times between the CG and cooling 
estimates of v^. 

At asymptotically large cooling times we expect N and 
V2N of the cooled configuration to vanish. If the CG val- 
ues of these do not vanish, then they are artifacts of the 
CG. We subtract the residual values from the correspond- 
ing values before cooling. The result is referred to as the 
corrected CG values. The cooling and the CG computa- 
tions are expected to agree at late times, as the system 
becomes increasingly weakly coupled. The two methods 
agree for N at fairly early times. For V2, this conver- 
gence occurs at much later times, since, as we shall see 
in the following, vi is dominated by very soft modes with 
momenta j?t < A s o. 

In Figure |^ we plot V2 reconstructed from the cooling 
time history of only the potential terms in T xx — T yy , 
along with the CG values (also including potential terms 
only) as a function of n c h/n to t for different values of A s0 i? 
as discussed in the figure. The systematic errors repre- 
sented by the band (for A s0 = 18.5 [p7l) are primar- 
ily due to limited resources available to study the slow 
convergence of the cooling and CG computations. We 
have studied the late time behavior of V2 for one impact 
parameter-the results are shown in the figure. 

The asymptotic values of V2, as predicted by the model, 
undershoot the data. This disagreement notwithstand- 
ing, our results show that a significant V2 can be gener- 
ated by the classical fields. For very peripheral collisions, 
where the gluon density may be too low to justify the 
classical approximation, the predictions of the model are 
not reliable. Interestingly, the dependence of v 2 on A s0 i? 
is rather weak. For a fixed impact parameter, the model 
predicts that, as A s rj-R — > 00, the classical contribution 
to the elliptic flow goes to zero. This is because increas- 
ing A s qR is equivalent to increasing R for fixed A s rj and 
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FIG. 2: V2 is plotted (for impact parameter b/2R = 0.5) 
versus A s r. Solid squares are cooling results using only the 
potential contribution. Solid diamonds include the contribu- 
tion of both potential and kinetic terms. The Coulomb gauge 
result, including both the potential and the kinetic contribu- 
tions, is shown in open circles. Open squares are corrected 
CG results-see text for an explanation. 
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FIG. 3: The centrality dependence of V2 at the earliest times 
in Fig. |^ is computed using cooling (open symbols) and CG 
(filled symbols). Results are for A s qR spanning the RHIC- 
LHC range, specifically, A s rj-R = 18.5 (squares), 37 (triangles), 
and 74 (stars). Full circles denote preliminary STAR data f7J. 
The band denotes the estimated value of V2 when extrapolated 
to very late times. "Corrected values" denotes the late time 
cooling and CG result for A 3 qR = 18.5 at one centrality value. 



therefore reducing the initial anisotropy. 

In Fig |, v 2 (pt) is plotted for b/2R = 0.75 for 
A s rji? = 74. Our calculations show that the elliptic 
flow rises rapidly and is peaked for p^ ~ A s0 /4 before 
falling rapidly. The theoretical prediction J2^| is that for 
Pt ^> A s0 , V2(pt) ~ A^o/Pt- The lattice numerical data 
appear to confirm this result-better statistics are required 
to determine the large momentum behavior accurately. A 
couple of comments about our result are in order. Firstly, 
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FIG. 4: «2(pt) as a function of transverse momentum in di- 
mensionless units for A s oR = 74. 



even though A^ is large, it may differ considerably from 
the color charge squared in the region where the nuclei 
overlap. This may explain in part why the momenta 
are peaked at smaller values of pt- Secondly, the domi- 
nant contribution of very soft modes to v 2 helps explain 
why the cooling and CG computations differ until very 
late times. The soft gluon modes have large magnitudes 
and therefore continue to interact strongly until very late 
proper times. Concomitantly, the occupation number of 
these modes is not small and the classical approach may 
be adequate to describe these modes even at the late 
times considered. 



DISCUSSION 

We now turn to the theoretical interpretation of the 
RHIC v 2 data in the CGC approach. It is clear from 
Fig. H that our result for v 2 contributes only about 50% 
of the measured v 2 for various centralities. Our px dis- 
tributions also clearly disagree with experiment || p9[ . 
Naively, one could argue that the classical Yang-Mills 
approach is only applicable at early times so additional 
contributions to v 2 will arise from later stages of the col- 
lision. While there is merit in this statement, it is also 
problematic as we will discuss below. 

The reason the situation is complex is as follows. We 
observed that it takes a long time r ~ R to obtain a 
significant elliptic flow. At these late times, one would 
expect that the classical approach would be inapplicable 
due to the rapid expansion of the system. On the other 
hand, we have seen that v 2 in the classical approach is 
dominated by soft modes which are strongly interacting 
and don't linearize even at time scales r ~ R. Clearly, 
the soft modes cannot be treated as on-shell partons even 
at times r ~ Rl This is the message one obtains from 
Fig | 



The correct way to treat the theoretical problem may 
be as follows. Hard modes with k t > A s linearize on very 
short time scales r ~ 1/Q S - Their subsequent evolution 
is treated incorrectly in the classical approach, which has 
them free streaming in the transverse plane. In actuality, 
they are scattering off each other via elastic gg — ► gg and 
inelastic gg <-> ggg collisions which drive them towards an 



isotropic distribution 30 . This dynamics would indeed 
provide an additional pre- equilibrium contribution to v 2 
and is calculable. An effect to consider here would be 
the possible screening of infrared divergences in the hard 
scattering by the time dependent classical field. More 
complicated is the effect of these hard modes on the clas- 
sical dynamics of the soft modes and on their possible 
modification of the contribution of the latter to v 2 ■ One 
has here a little explored dynamical analog to the in- 
terplay of hard particle and soft classical modes in the 
kinetic theory of Hard Thermal Loops |n| . 

Were the system to thermalize, both these effects 
would complement the hydrodynamic component of el- 
liptic flow [pi- A quantitative study of the anisotropy 
generated in this intermediate regime would therefore be 
very useful in our theoretical understanding of the data. 

Finally, we note that v 2 is extracted only indirectly 
from a variety of techniques-in particular, two and four 
particle cumulant analyses |52|. Recently, it has been 
proposed that non-flow two particle correlations explain 
the v 2 data |33). It is unclear whether this model 
can explain other features of the measured azimuthal 
anisotropy. In our approach, a procedure very similar 
to the experimental approach can be followed and two 
and four particle correlations can be determined. This 
work will be reported in the near future. 
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